Closed, spirograph-like orbits in power law central potentials 
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Bertrand's theorem proves that inverse square and Hooke's law-type central forces are the only 
ones for which all bounded orbits are closed. Similar analysis was used to show that for other central 
force laws there exist closed orbits for a discrete set of angular momentum and energy values. These 
orbits can in general be characterized as "spirograph" -like, although specific orbits look more "star"- 
like or "triangular." We use the results of a perturbative version of Bertrand's theorem to predict 
which values of angular momentum and energy result in closed orbits, and what their shapes will 
be. 
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be found at http://scitation.aip.org/ajp/. 



I. INTRODUCTION 



O |. It is well known that there are only two central force laws for which all bounded orbits are closed. ^'^ By "closed," 
■ we mean that the orbiting object returns to the same spatial location with the same velocity in a finite amount of time 
, (specifically, it returns to the same location in phase space). By "bounded," we mean that the distance r between 
the orbiting object and the central object always remains between two fixed values, rmin < J" < I'max, called the 
^ \ radial turning points, or in the case of elliptical orbits they are called periapsis and apoapsis. This result is known as 

. Bertrand's theorem,— first obtained in 1873. 
O ' If the attractive force is represented by a power law, F = —kr^, then only A — —2 (an inverse square force given 
' ^ , by Newton's law of gravitation or Coulomb's law) and A = 1 (a spring-like force given by Hooke's law) admit closed 
orbits, both of which happen to be elliptical. In fact, the orbits in these two potentials satisfy the additional criterion 
that they are "non-crossing." A bounded orbit does not cross itself in configuration space if the ratio of its orbital 
period, Tort, to the period of its radial oscillations. Trad, is an integer. Here, we use the parameter /3 = Torb/Trad 
' to denote this ratio. In the case of Newtonian gravity /? = 1, which means that there is only one periapsis and 
y—^ • one apoapsis per orbit, and the central body resides at one focus of the elliptical orbit. In the case of Hooke's law, 
K*" ' however, (3 — 2, and the central body is located at the center of the elliptical orbit. This means that there are four 
0^ turning points (two close, two far) in each orbit. 

The implications of Bertrand's theorem have been investigated extensively, ranging from the symmetries inherent 
in the potentials to the deep connections between classical and quantum mechanics that it reveals.— The fact that an 
orbit is closed means that, besides energy and angular momentum, there must be an additional conserved quantity — 
25 the Runge-Lenz vector<^ Also, closely related to the fact that only A — —2 and A = 1 admit closed classical orbits is the 
result that these two potentials result in an exactly solvable Schrodinger equation;^ In addition, these two potentials 
are "dual" in the sense that one problem can be obtained from the other by a change of variable)^ In addition, many 
authors have obtained proofs of the theorem that are more elegant and pedagogical than the originalp'"— and central 
potentials other than power-law have been investigated.— 

In this work we focus on analytical methods suitable for the intermediate mechanics student, as well as numerical 
5J] ' techniques that can be used to find closed orbits (especially those with interesting shapes) in central forces other than 
inverse square or Hooke's law. The types of closed orbits that can be obtained are introduced in Sec. [Hi In Sec. IIII[ 
Bertrand's proof of his eponymous theorem is briefiy outlined, and a more pedagogical proof, first given by Brown;^ 
is covered in detail. This detail is needed because Brown's method includes the mathematical insight necessary 
to analyze large amplitude perturbations from stable circular orbits. Finally, in Sec. IIVI we obtain conditions that 
must be satisfied so that these large-amplitude orbits are closed, and several representative trajectories are obtained 
numerically. 



II. CLOSED ORBIT CLASSIFICATION 



For all power law central forces other than inverse-square and Hooke's law, most orbits, while they may remain 
bounded, are not closed. However, there are three cases in which the orbits are closed. First, as long as A > —3, all 
power law central forces exhibit a stable, closed, circular orbit at the radial location ro where the effective force F^g 
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Feffiro) = -kr^ + — = 0, (1) 
mrg 

where the second term is the "centrifugal force," L is the (constant) angular momentuni, and m is the mass of the 
orbiting object. This stable radial location is given by = L'^/mk. Of course, if A < —3 there are no stable 

circular orbits. Indeed, when A — —3, Newton showed that the trajectory is a so-called Cotes' spiraljiS, 

Any deviation from a circular trajectory allows the possibility that the orbit may no longer be closed. However, the 
second case in which closed orbits arise is when the orbiting object is perturbed only slightly from the stable circular 
orbit. If r is infinitesimally close to rg, then the effective potential energy (defined as usual by F^ff = ^dUeff/dr) can 
be expanded about rp 

Ueff^U{ro) + \k,ff{r-rof, (2) 



where the effective spring constant is 



fc,jy=(A + 3)-^. (3) 

Two types of periodic motion are now superposed, the previous orbital motion as well as a radial oscillation in the 
simple-harmonic effective potential of Eq. ([2]). The period of these radial oscillations is Trad = '^T^\fnijk^. Since the 

orbital period, obtained by taking a ratio of the circumference, 2T:ro, to the orbital velocity, ^Jlzr^^^^^~fm, is given by 

/ i-> l-~\-P ^-Vl/^ ^-TTT/~l l-li-M-" 1 /~1 /-I CI 1 ClU 



Torb = 27r-i/ mrQ /k, the ratio of the two periods i 



/3^±i^ = x/AT3. (4) 

^ rad 

If A is such that Torb/Trad is a rational fraction, (3 — p/q, where p and q are integers, then this "almost-circular" orbit 
will be closed. Hence, for certain forms of the power law that satisfy 

A=%-3, (5) 
9 

orbits that are only slightly perturbed from a circular orbit are closed. Bcrtrand^ was able to show that in two special 
cases, p — q and p — 2q, corresponding to A = —2 and A = 1, respectively, orbits with large (not just infinitesimal) 
deviations from a circular trajectory remain closed. 

This analysis suggests that other solutions of Eq. ([S]), e.g., A = 6 {p — 3q) and A = 13 (p = Aq), admit closed orbits 
for infinitesimal perturbations from a circular orbit. However, we show in Sec. llVl that for all values of A that exhibit 
stable circular orbits, finite perturbations from a circular orbit can result in values of /? that are rational fractions. 
This is the third case, mentioned above, in which the orbits are closed. In this situation, most of the energy-angular 
momentum parameter space results in non closed orbits, but certain discrete values of these two parameters result 
in closed orbits. Further, most of these are "crossing" orbits in which the trajectory crosses itself one or more times 
before returning to the original location, which means that they correspond to a rational fraction where q 1. There 
are a few instances in which (3 is an integer, though, and these orbits can be triangular (/3 = 3) or even square (/? = 4) 
in shape. These large values of /3, however, require large positive values of A. 

A. Central forces that are not power laws 

Of course, the previous conclusions, as well as the analysis below, are not limited to power law central forces. 
Gauss's Law implies that an arbitrary (but spherically symmetric) mass density distribution p{r) results in a central 
force law 

Fr-^~\ f p{r'y^dr'. (6) 
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For any particular density distribution of interest, the effective force and potential energy, the stable circular orbit 
radius, and the ratio of the orbital and radial oscillation periods, Eqs. can all be obtained. Any parameters 

describing p(r) will of course replace A. For power law central forces, the self-consistent density distribution is 

p[r)^{\ + 2y-\ (7) 

In fact, quite a bit of theoretical work has been done on the problem of orbits in the gravitational potentials of galaxies 
and globular clusters.— For example, Adams and Blocbi^ analyzed orbits in the so-called Hernquist potential 

^ (8) 



(r + rs)2 

where is the length scale of the potential, and the potential is due to an extended mass distribution with density 
p(r) ~ l/r(r -|- Ts)^. This distribution turns out to be a good approximation for elliptical galaxies and dark matter 
haloes. The focus in these studies has been on understanding how the orbits affect the dynamics of the system, and 
not on whether each individual orbit is closed or not. Also, Strucki^ was able to analytically solve for the orbits 
using the so-called "epicycloid" approximation, which assumes the orbit is a precessing ellipse whose shape can be 
expressed as a function of the type 

J{e)[l + ecos{l-b)e] (9) 

where e is the eccentricity and b determines the precession rate. Of course, the parameters e and 6, along with the 
function f{9) must be determined from the form of the potential. This technique allowed him to obtain the result 
in Eq. ^ above, and therefore obtain orbital resonance conditions that can assist understanding galactic dynamics, 
such as bars in spiral galaxies. 

III. FINITE RADIAL PERTURBATIONS 

Bertrandi used the well-known orbit equationsi^ to express AS as an integral over the radial excursion r, where 
AS is the angle swept out by the trajectory. In order for the orbit to be closed, he then required that this integral, 
when evaluated between two neighboring turning points, be a rational fraction times tt, or in our notation, 

A^-^=P , (10) 

where r„jj„ and r^ax are roots of the denominator. He took a global approach, simultaneously expanding the integral 
for small oscillations about a stable circular orbit as well as letting r^ax ~^ oo. He was then able to show that the 
requirement in Eq. (jlOl) means that V{r) must be a power law with A = — 2orA=l. Unfortunately, his proof does not 
easily divulge any physical insight. On the other hand. Brown's method,- in which he solved for the periodic motion 
in the anharmonic potential [see Eq. (jl3p ] near the radius of the stable circular orbit, not only proves Bcrtrand's 
theorem, but also allows the derivation of a closed orbit criterion that is valid for any power A. 

Here we outline Brown's method, and quote the results that are relevant to the present discussion. First, he solved 
the dynamical equation for radial motion in the potential given by Eq. ^ by assuming that the object is in an initially 

stable, circular orbit with r — vq and orbital speed vg^ — \J kr^^^ /m. Then a small radial impulse is imparted to 

the object (in order to conserve the angular momentum L) which results in a nonzero radial velocity Wr,o- Of course, 
the subsequent trajectory consists of a harmonic oscillation of the radial coordinate, r, 

r{t) ^ ro + e cos LOot, (11) 

where Wq = keff/m is just the frequency of small radial oscillations, as we obtained above, and e is the amplitude of 
the radial oscillations. There is a simple relation between the initial radial velocity w^.o a-nd the amplitude e, which is 
Vr,o = Woe, or 



^ = yAT3 -), (12) 



and which comes from elementary simple-harmonic-motion theoryJi 
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Bertrand's theorem, however, is a statement about the character of finite radial oscillations, and the restriction to 
infinitesimal amplitudes must therefore be relaxed. It turns out that it is sufficient to retain one more term, the cubic 
term, in the expansion in Eq. ([2]), which becomes 



U. 



eff' 



U{ro) + —xir- r^f + - 



12) 
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(13) 



and then apply the classic solution to this anharmonic oscillator problem, which was given by Landau and LifshitzJ^ 
The technique consists of seeking a solution that is a series of "successive approximations." The first order approx- 
imation is just Eq. (jlip . while the second and third order approximations include oscillations at harmonics of the 
fundamental frequency, cos 2ciji and cos Suit. Here, w is the exact anharmonic oscillation frequency, slightly shifted 
from loq by a term that is proportional to 



(A-10)(A-1) 
12 



ro 



(14) 



A well-known example of this effect is the large-amplitude pendulum, whose exact restoring force is proportional to 
sin 9 ^6 — 9^/31, and an inclusion of the cubic term results in an amplitude-dependent period. 

Brown also showed that for large amplitudes the angular velocity of the orbital motion is also slightly shifted by a 
term that is proportional to e^. Using his notation 



(0)' - Wo 



1 + (A-1) 



ro 



(15) 



where 9 is the angular orbital velocity, the brackets () indicate an average over one orbital period, and the subscript 
denotes the stable circular orbit value in the limit e — > 0. Retaining only terms of lowest order in e, the ratio of the 
two periods is 



rp2 

orb 
rad 



{9Y 



(A + 3) 



(A-l)(A + 2) 
12 



(16) 



This is Brown's main result, and it proves Bertrand's theorem. For all orbits to be closed, the ratio of the two periods, 
/3, must be independent of the radial amplitude, and this is only true when the coefficient of is zero. That is, A = 1 
or A = —2, as previously stated. As it must, Eq. (|16|) also contains the limit given in Eq. (|4]), which might be called 
a "restricted version" of Bertrand's theorem: "For infinitesimal perturbations, e — )• 0, Eq. (IT51) reduces to Eq. (jH), 
which means that the condition for closed orbits is Eq. ([S])." 

From a practical perspective, however, to integrate Newton's second law numerically and obtain a trajectory, it 
is the initial conditions, v^^q and ve.Oi that must be specified. In addition, it is the parameter (3 that is of primary 
interest, not the radial amplitude e. It is useful, therefore, to eliminate e from Eq. using Eq. ([T^. to obtain 



« (A + 3) - 



(A - l)(A + 2) (Vr,0 
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(17) 



Trajectories that demonstrate the restricted version of Bertrand's theorem (obtained by numerically integrating New- 
ton's second law using a Runge-Kutta 4th order method) are shown in Figs. [1] and [5] for a force law parameter A = 6. — 
In Fig. [U an initial condition of Vrfl/vgg = 0.003 results in an almost circular orbit. Since the initial radial velocity 
is small, the radial amplitude is likewise small, and Eq. (IT^ predicts e/ro « 10~^ for the parameters chosen, which 
agrees with the numerical result shown in Fig. [TJb). In addition, Eq. ([4]) predicts /3 = 3, which is also seen clearly 
in Fig.fljb), even though the radial oscillation is not perceptible in Fig.[lja). All trajectories in this paper share the 

following initial conditions: = tq, = 0, and Vy Q = vg Q = \J kr^^^ /m. This means that if Vr,o = 0, then the orbit 
is stable and circular. It also means that the angular momentum remains fixed. Varying the initial radial velocity 
changes the orbit shape because the total energy varies. 

When the radial impulse imparts a large radial velocity, say Vrfl/vgfl = 0.3, the closed nature of the orbit is lost, 
even though it is still bounded. This can be seen in Fig. [51 For A = 6, Eq. becomes 
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Vr,0 
V0fi 



(18) 
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or /3 « 3— (5/9) (wr,o/''^e.o) , which shows that the orbital period decreases to less than three times the radial oscillation 
period as the radial amplitude increases. This is indicated in Fig. [^b) by the fact that the radial position does not 
quite return to tq after one complete orbit. We can confirm this mismatch quantitatively using Eq. (|f 8|) . which gives 
/3 « 2.95, and this means that when 9 — 27r, the radial oscillation should have a phase 27r x 2.95, and a displacement 
of 0.1sin(27r x 2.95) — —0.032, and this is just what is observed in Fig. [2ljb). The amplitude is also consistent, for 
Eq. (fT2|) predicts e/rg « 0.1, which again agrees with the numerical result in Fig.[2][b). 



IV. CONDITIONS FOR CLOSED ORBITS 



Now that we have confirmed numerically the restricted version of Bertrand's theorem, along with the fact that the 
orbit does not remain closed when the radial amplitude is not infinitesimal, we can now investigate the conditions 
that allow large amplitude orbits (in power laws other than A = — 1, 2) to be closed. In fact, Eq. ([T7|) is just such 
a condition. Above, we used Eq. (|17p to predict the value of f3 (and whether it is a rational fraction or not) from a 
knowledge of the initial conditions (e.g., Wr,o) and it worked as long as fr,o was small. Now, however, it is clear that 
Eq. pT| also indicates that there can be closed orbits for any value of A, as long as Vrfl has the correct value. To see 
this, invert Eq. ([TT)) to obtain Vrfi as a function of /3 

\ 2 

'''•■"^ - (A + 3-/32). (19) 



ve,oJ (A-l)(A + 2) 

In this case, we first choose the force law parameter A and then the desired ratio of the periods, f3. Then, Eq. (fT9|) 
predicts the initial radial velocity needed to obtain that particular closed orbit. Of course, the larger that the difference 
is between /3 and VA + 3, the larger the radial oscillation, and Eq. (jl9p represents a poorer approximation. 
For example, again considering the force law parameter A = 6, Eq. (jl9p reduces to 

-.o^^ 3(9-/3^) ^2^^ 



V0fij 10 

It is clear that /3 = 3 is the small radial oscillation limit since it predicts an initial radial velocity of u,.^o = 0. In 
addition, since Wr,o raust be real, (3 will always be less than 3 as the orbit deviates from a stable circle. This fact was 
already clear from Eq. ()16|) . As Vr,o increases from zero, f3 will take on a continuum of real values less than 3, most 
of which will not be rational. However, /3 will pass through an infinite number of discrete values that are rational, 
implying that the corresponding orbit will be closed. For the trajectory in Fig. [21 /3 ~ 2.95, and it is probably not 
rational, since it was obtained by fixing Vr^. A rational fraction near this value is /3 = — which means that 
the orbit will have 59 radial oscillations for every 20 orbits about the center. Such an orbit is shown in Fig. |31 which 
is clearly closed with the correct value of /?. However, the initial radial velocity needed to obtain this orbit is not 
quite the prediction of Eq. ([20)) . which is Vrfi/vgfi ~ 0.29875. This is because although Eq. ((20)) follows from Eq. p6t . 
which is valid for large enough radial amplitudes to prove Bertrand's theorem, it represents a poorer approximation 
as (3 — P) increases. In order to determine the correct value of Vr,o needed for such an orbit, a more sophisticated 
numerical technique is required. 



A. Numerical determination of closed orbits 



There are two methods that can be used to find the necessary value of Vr.o that results in an orbit with a particular 
/3: brute force search and root finding. Both methods can successfully utilize the technique of Poincare's surface of 
section)^ which takes the continuous time evolution of a high-dimensional trajectory and replaces it with a discrete 
mapping in fewer dimensions, usually two. In the present case, we plot in radial phase space (i.e., f versus r) the 
locations where a particular trajectory crosses the positive x axis, for example. Then, closed orbits can be found when 
the trajectory returns to the same phase space location after an integral number of orbits. The surface of section for 
the trajectory in Fig. |3]is shown in Fig.U] Since the initial conditions were r = ro and r = Vr.o, where Vr^o is positive, 
the initial location in Fig. 2] is denoted by a circle. After 20 orbits, and therefore 20 crossings of the positive x axis, 
the trajectory returns to the same phase space location. This confirms that the orbit is closed. In fact, the trajectory 
can be followed for several "recurrence periods," i.e., 40 or 60 orbits, to make sure that the Poincare section is periodic 
in the long term. 

In addition to confirming that the orbit is closed, the surface of section suggests a technique that works for the 
second method: root finding via the shooting methodi^ Here, the shooting method works in the standard way, by 
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casting the problem as a two-point boundary value problem. The initial condition is varied — in this case Vr.o (the 
other three initial conditions, Xq, uq, and vg^, remain fixed) — and the equation of motion is integrated until the 
desired final condition is obtained. The final condition here is that for an orbit with /3 — ^, the distance in phase 
space between the "zeroth" crossing of the positive x axis and the qth crossing be zero, i.e., they must be identical. 
Of course, a good initial guess for Vrfi is needed, and this is supplied by Eq. In addition, a robust root- finding 

method must be employed. Since the derivative of our "function" (distance in phase space as a function of Wr,o) 
is not available analytically, and since the tolerance of the root-finding method should not exceed the tolerance of 
the numerical integration, the simple secant method should work fine. On the other hand, since the phase space 
distance is a positive definite quantity, the desired distance is not just a root, but also a minimum. For this reason, a 
minimization method, such as Brent's method, can also be used. It turns out that in practice, either method works 
fine. 

In principle, orbits with any allowed value of /3 can be found provided the initial guess for Vr^o is accurate enough. 
Even if Eq. does not supply a sufficiently accurate first guess, the "distance function" versus Vr,o can easily 
be calculated and plotted, and a better first guess obtained. For A = 6, several closed orbits were found using this 
method, and the values of f3 and Ur,o for each orbit are shown in Fig. \E\ The small amplitude relationship, Eq. (|20|) . 
is also shown, and it can be seen that the two deviate when the radial amplitude becomes large. 

What do these large amplitude orbits look like? Besides the orbits with large values of q, which are close to 
circular, the crosses in Fig. [5] indicate a few orbits with small values of q (of course with (3 still less than 3). The 
orbit with the smallest value of q is /3 = |. The initial radial velocity and amplitude predicted by Eqs. (PH)) and 
(IT^ are Wr.o/^e.o = 0.908, and e/j'o = 0.303. Since this radial oscillation amplitude is large, the small amplitude 
result in Eq. (jl6p is not applicable. A search of parameter space (using the secant method explained above) reveals 
that the necessary initial radial velocity is v^fi/vg^o — 1.36671, and this orbit is shown in Fig. [51 Even though it can 
be classified as spirograph- like, because /? is a ratio of two small integers (and is greater than unity), the orbit has 
the appearance of being more "star"-like. Other similar values of /?, for example the cross labeled | in Fig. O are 
consistent with orbits that also have a star-like appearance. The radial displacement of the star-like orbit turns out 
not to be centered on rg, which is to be expected from a large-amplitude, anharmonic oscillator. A rough estimate 
from the numerical solution gives e/ro « 0.405, which is significantly larger than that predicted by Eq. (|12|) . 

B. Large amplitude orbits and non-crossing orbits 

For a given value of A, what is the range of possible values of /3? We have seen that for A = 6, /3 must remain less 
than three. For other values of A, /? is also restricted, and this restriction is determined by Eq. (|17|) . which shows 
that /? must be either greater than or less than -y/A -I- 3 depending on the sign of the coefficient oi v^q. For A = 6, the 
coefficient is negative, — (A — l)(A-|-2) = —40 < 0, which means that /3 < 3, as we have already discovered. This result 
divides the parameter space into three regimes, and the boundaries between these regimes are just the two special 
cases of Bertrand's theorem: 

regime I A < -2 /3 < v^ATS 

regime II -2 < A < 1 /3 > ^AT3 (21) 
regime III 1< A /3 < ^/AT3 

Regime I is actually restricted to —3 < A < —2 because we are only interested in bounded orbits. The character of 
these orbits can be different from those we have already studied, because /3 will always be less than one. In fact, as 
we have defined it, /? must be positive definite, so for regime I it must lie in the range < (3 < VA -I- 3. If A = —2.5, 
for example, (3 < l/-\/2 « 0.707. In this regime, when /? is a ratio of two fairly large integers, the orbits are similar to 
the orbit in Fig. [31 To see this, the case of A ~ —2.5 and (3 = is shown in Fig. [T] The only difference in character 
between the two orbits is that in Fig. [7] the number of orbits is greater than the number of "furthest approaches" 
where r = Vmaxi rather than vice- versa. On the other hand, regime I allows a new type of orbit because when (3 is the 
ratio of two small integers, the nature of the trajectory radically changes. Again, for A = —2.5, the closed orbit where 
/? = I is shown in Fig. |51 Because only two values of rmin (and rmax) can occur during the course of three orbits, 
the trajectory looks very different from Fig. [3] In fact, this orbit appears more "loop"-like than spirograph-like. This 
character comes from the fact that /3 is the ratio of two small integers and is less than one. Regime I is the only case 
where /3 can be less than unity. 

In addition to non-power law forces, Strucki^ focused on power-laws in regime II, because these describe galactic 
potentials well. He showed that in addition to the criterion in Eq. (|2ip . (3 was restricted to /3 < 2. The global analysis 
of Bertrand^ also reveals this fact, and in particular shows that /3 — >■ 2 in the limit that rmax ^ oo. (For the numerical 
solutions in this study, this limit corresponds to I'r.o/^^e.o ~^ oo.) This limit also explains why the exact numerical 
solutions in Fig. [S] are all in the range 2 < /3 < 3. The allowed values of (3 for all three regimes are shown in Fig. [HI 



7 



Because of the restriction 1 < /? < 2 in regime II, the types of orbits have the same character as Fig. |3l That is, they 
are of the spirograph type, and except for A — —2 and A = 1, they cannot be non-crossing. 

The final type of orbit with an interesting character occurs only in regime III. These are characterized hy q — 1, 
which means that the orbit is non-crossing. The case of A = 6, studied above, does not admit a large-amplitude, 
non-crossing orbit, since P = j is restricted to the nearly circular case, and (3 — j is not accessible with a finite value 
of Vr.o- However, if VA + 3 > 3, then a finite amplitude orbit can be consistent with 13 ^ j, resulting in a closed, 
non-crossing, "triangular" -shaped orbit. This is shown in Fig. [TU] for the force law parameter A = 8. Any value of A 
greater than 6 will, of course, admit a triangular orbit if Vr.o has the proper value. "Square" -shaped orbits can also 
occur when VA + 3 > 4, i.e., A > 13, and one is shown in Fig. [TT] where A = 20. Higher order "polygonal" orbits are 
also possible, but they require increasingly larger minimum values of A. 



V. CONCLUSION 



Closed orbit trajectories of several different types have been found for central forces that are of a power-law type. 
Besides the well-known elliptical orbits that arise from Coulomb's law r^^) and Hooke's law r), we have shown 
that closed orbits exist for all power law central forces, Fr — ~kr^ , when A > —3. Over the largest part of parameter 
space, the closed orbits are spirograph-like, with many self crossings before they return to their original location. 
However, when ^ is a ratio of two small integers, then the orbits become more "star" -like (Fig. [6]) or "loop" -like 
(Fig. [5]). Finally, non-crossing orbits, when /3 is an integer, occur for large values of A, and can be triangular, square, 
or polygonal. 
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FIG. 1: (a) Trajectory in the orbital plane (x, y) for the case A = 6, with the initial conditions .xo = ro, yo = 0. and 
Vr,o/ve,o = 0.003. (b) Radial location, r = y/x^ + y^, as a function of angular position 9 = arctan(t//x). The three radial 
excursions of the first orbit can be clearly seen, with an amplitude e/ro = 10~^, as well as the fact that when = 27r, r returns 
to ro. 
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FIG. 2: Same as Fig. [T] except that Vrfi/vefi = 0.3. This implies that e/ro ~ 0.1, and j3 ~ 2.95. This means that the period of 
radial oscillations is slightly more than one third of the orbital period, so that when 6 = 2n, the radial oscillation should have 
a phase 2it x 2.95, and a displacement of 0.1 sin(27r x 2.95) = —0.032, as observed. 
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FIG. 3: Same as Fig. EJa), except that Vr,o/vB,o = 0.30748. The orbit is closed, with /? = ||. The value that the small 
amplitude approximation, Eq. (|20|l . predicts for the radial impulse is Vr,o/vg,o ~ 0.29875. However, the radial amplitude is too 
large, e/ro ~ 0.1, for this prediction to be exact. A numerical root finding search was used to obtain the correct value. 
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FIG. 4: Poincare surface-of-section for the parameters in Fig. [S] The trajectory's location in {r,Vr) phase space is plotted 
each time it crosses the positive a:-axis. Forty revolutions about the center are shown, so that each cross is really two crosses, 
depicting subsequent passings. The fact that the crossing locations are identical means that the orbit is closed. 



13 




FIG. 5: Relationship between /3 and Vr,o for A = 6 from Eq. (|20|l (solid line) and exact numerical calculation for closed orbits 
(crosses). The values of /3 corresponding to the five largest amplitude closed orbits are also indicated. 
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FIG. 6: Same as Fig. O except that «r,o/«fl,o = 1.36671. The orbit is closed, with /3 = |, and the radial amplitude is 
approximately e/ro ~ 0.405. 
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FIG. 7: Same as Fig. O except that A — —2.5, Vr,o/v0,a = 0.25034. The orbit is closed, with (3 = j^. Again, the anharmonicity 
of the potential means that the radial amplitude is not symmetric about ro. The numerical result gives r-min ~ 0.743ro and 
r-max ~ 1.581ro. 




FIG. 8: Same as Fig. [T] except that Vr,o/ve.o = 0.48733. The orbit is closed, with l3 = ^. 
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FIG. 9: Allowed values of /3 for the three regimes. The curve is the linear relationship between /3 and A for nearly circular 
orbits, given by Eq. Q. The shaded regions indicate allowed values of (3 for large radial oscillations in the three regimes, whose 
limits are given by Eq. (|2ip and in the subsequent text. 




FIG. 10: Same as Fig. El except that A = 8 and Vr,o/ve,o = 0.68492. The orbit is closed, with P = \. 




FIG. 11: Same as Fig.[TOl except that A = 20 and Vr,a/vgs) = 0.54372. The orbit is closed, with 13 = \. 



